Non-universal results induced by diversity distribution in coupled excitable systems 
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We consider a system of globally coupled active rotators near the excitable regime. The system 
displays a transition to a state of collective firing induced by disorder. We show that this transition 
is found generically for any diversity distribution with well defined moments. Singularly, for the 
Lorentzian distribution (widely used in Kuramoto-like systems) the transition is not present. This 
warns about the use of Lorentzian distributions to understand the generic properties of coupled 
oscillators. 
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Synchronization phenomena play a prominent role in 
many branches of science [l|. They have been analyzed 
in terms of phase models which successfully describe sys- 
tems of weakly coupled limit cycle oscillators. In partic- 
ular, the Kuramoto model [2[ has become a paradigm for 
the study of synchronization (for reviews see [l|,|3|,|4l|). It 
shows how synchronized behavior can appear when the 
competitive effects of coupling and diversity among the 
individual units are present. 

Amongst other possible sources, diversity in the oscil- 
lators is usually introduced by taking their natural fre- 
quencies from a probability distribution. Although, on 
general grounds (central limit theorem), this distribution 
should be well approximated by a Gaussian form, theo- 
retical studies usually consider a Lorentzian form since 
it allows for an easier analytical treatment. It is gener- 
ally believed that the main results concerning the global 
synchronization properties are qualitatively independent 
of the precise form of the distribution as long as it is 
symmetric and unimodal. In this Letter, however, we 
show that a variant of the Kuramoto model displays or 
not a reentrant diversity-induced transition into a state 
of collective firing, depending on the type of distribution 
used. This transition is present in all the (symmetric 
and unimodal) distributions studied, but not in the case 
of the Lorentzian. The non-generic behavior of the sys- 
tem with a Lorentzian distribution of natural frequencies 
warns about the indiscriminate use of some recently pro- 
posed methods [5[ in order to understand generic prop- 
erties of coupled oscillators. 

We consider the following variant of the Kuramoto 
model which describes the dynamics of an ensemble of 
globally coupled active rotators <f>j(t),j = 1, ..., N [6|: 
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density function g(uj), with mean value ZJ and variance 
a 2 . The model is equivalent to the Kuramoto model with 
zero average frequency and an external periodic driving 
of frequency — 57, as it can be easily seen with the change 
of variables <j)j — > </)j — ZJt. Throughout the paper, be- 
sides the well-known Gaussian and uniform distributions, 
we will be considering a general family of distributions 
L™(uj), for n > 0,mn > 1, defined as: 
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The variance is finite only for mn > 3 and it is given 

The Lorentzian distribution 
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by a 1 = A 

corresponds to L\(uS) and has, hence, an infinite variance, 
although we still will use A as a measure of diversity. 

To characterize the collective behavior we use the time- 
dependent complex variable r(t) [2j,|6|: 
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The Kuramoto order parameter p = (\r(t)\), where (• • • ) 
denotes time average, is known to be a good measure 
of collective synchronization in coupled oscillators sys- 
tems, i.e. p ~ 1 when the oscillators synchronize ((f>j ~ 
<pl,Vj,l), and p ~ for desynchronized behavior. 

For To < 1 the system displays three different regimes: 
(i) for small diversity, almost all units are at rest at simi- 
lar fixed points (r(£) does not depend on time and p~l); 
(ii) increasing diversity one enters a dynamical state in 
which a macroscopic fraction of units fire at (roughly) 
the same time (r(i) is time dependent while p is still close 
to 1); (iii) for even larger diversity, the system reenters 
a desynchronized state (p small and r(t) time indepen- 
dent). To discriminate between static entrainment and 
collective firing, regimes (i) and (ii), we use the order 
parameter introduced by Shinomoto and Kuramoto [7|: 



A natural frequency ujj < 1 (respectively, u)j > 1) cor- 
responds to an excitable (respectively, oscillatory) be- 
havior of the rotator j when it is uncoupled. K is the 
coupling intensity. Diversity is introduced by considering 
that the u/j's are distributed according to a probability which differs from zero only for collective firing. 



C = <|r(t)-(r(t))|), 
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An approximate theory to describe these three regimes 
was developed in [8] . The theory was independent of the 
form of the natural frequencies distribution and was also 
applicable to identical units subject to noise. A recent 
method developed by Ott and Antonsen (Jj |9( allows to 
solve exactly this model (and a large family of related 
ones) in the infinite number of oscillators limit and in a 
number of cases that include the Lorentzian g(u)). Childs 
and Strogatz [10| used this method to obtain the full bi- 
furcation diagram of the complex variable r(t) for the 
Lorentzian distribution. Contrarily to the results of [8j, 
their exact solution implies that there is no transition to 
collective firing increasing the diversity for ZJ < 1. The 
non-existence of the transition can be derived from the 
bifurcation diagram in the T3 — A space obtained using 
the ideas of [lOj, see FigfT] Regime (ii) takes place for the 
parameter region located to the right of the SNIC (saddle 
node on the invariant circle) bifurcation line and below 
the Hopf line. For ZU < 1 increasing A one never enters 
in this region, so there is no diversity-induced transition 
to collective firing. This situation is generic for all values 
of K, since it can be shown that the SNIC bifurcation 
always starts at W = 1, A = with positive slope. The 
model was also studied for the Lorentzian case with a dif- 
ferent approach in [11| and the same results were found. 
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FIG. 1: Bifurcation diagram of the model defined by Eqs. JT]) 
for a Lorentzian distribution of frequencies and K = 5. There 
is also a saddle-loop bifurcation line, not shown, that goes 
from the Takens-Bogdanov point (circle) to the saddle-node 
separatrix loop point (square) where SNIC line starts \W§ - 
Collective firing takes place in the region below the Hopf and 
right of the SNIC. 

We will give now the main sketches of the Ott- 
Antonsen method. Quite generally, we will show that, 
in fact, the method can be successfully applied to any 
non-singular distribution g(uS). Let /(oj, 0, t) be the den- 
sity of oscillators with frequency w and phase <f> at time 
t. This function obeys the continuity equation (conser- 
vation of the number of oscillators) : 
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of the Fourier expansion: 
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(c.c. denotes complex conjugate), satisfy the ansatz 

f n (w,t) = a(u,t) n , (7) 

then a(w,t) and r(t) — J dcj f d(j)e l ^f(u},4>,t) obey the 
integro-differential equations: 

-^-+iuja+-{[Kr + l}a 2 -Kr* - 1} = 0, (8) 



r(t) = / doj a(cu,t)* g(ui). 



with 4>(uj, 4>,r) — ui — sin <fi + 5ft(re ^). If the coefficients 
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The manifold defined by is invariant under the system 
evolution, so if that condition is fulfilled by the initial 
condition, it is fulfilled afterwards. Moreover, in [9j it is 
shown that for a Lorentzian g(u>) the long time evolution 
of r(t) is always described by this reduced manifold. 

If g(u) has a finite set of poles 0)1,0)2, • ■ • outside the 
real axis (as is the case for L™{uS) for even n and integer 
m, including the Lorentzian L\{lu), and a(u>,t) satisfies 
certain analyticity conditions, one can obtain ([9]) by con- 
tour integration. Then r(£) can be written is terms of 
ak(t) = a(u)k, t) and one can obtain a closed set of ordi- 
nary differential equations for au{t)- In the case of poles 
with multiplicity larger than one, r(t) depends also on the 
partial derivatives with respect to w, &%"(£) = a^ s '(uk,t). 
Equations for these new functions a s k (t) can be obtained 
by differentiating Eq.© with respect to u. For an ar- 
bitrary distribution g(cu), we can obtain an approximate 
evolution of the system by evaluating integral ([9]) using 
a finite, though large, set of values of uj and integrating 
numerically ([8]) for each one of these frequencies. 

In Figs. [5] and [3] we show the order parameters as 
a function of the diversity for several frequency distri- 
butions, obtained by direct simulation of Eqs. ([T]) and 
by applying Ott-Antonsen's method as described above. 
For the Gaussian and uniform distributions, but not for 
the Lorentzian, the regimen of collective firing (signaled 
by a nonzero value of the parameter £) is present for 
intermediate values of the diversity. Beyond the cases 
shown in the figures, the transition is also present for 
other symmetric distributions such as the exponential 

g ^ = 772 e ~^ lu ~ Wl ' or the family L ™ ( w ) for mn - 3 - 
Even £3(0;) which has infinite variance (but well-defined 
first moment) presents the reentrant diversity-induced 
transition (for values of a; close enough to one). Also, 
if we truncate the Lorentzian distribution at some finite 
value of uj, i.e. set g{uj) — 0, if \w — 5J| > C, the system 
shows this reentrant transition (we checked for C=50A). 
Furthermore, Fig. [3] shows that for finite-size systems 
with a Lorentzian distribution a maximum in £ is in- 
deed present, being quite visible up to a few thousand of 




FIG. 2: Order parameters p and ( as a function of the di- 
versity a coming from numerical simulations of the full sys- 
tem of Eqs.JT]) for N = 10 units in the cases of a Gaussian 
(dots) and uniform (crosses) distribution of frequencies. In 
each case, the solid lines are the result of the application of 
the Ott-Antonsen method using 10 4 values of w for the nu- 
merical integration of Eq.©. K = 5 and uJ = 0.97 



units. In fact, Lorentzian distributions in systems with a 
finite number of units arc effectively truncated, trunca- 
tion that disappears in the limit N — y oc. We conclude 
that the existence of the transition to common firing is a 
truly generic phenomenon and the results obtained using 
a Lorentzian distribution in the infinite system size limit 
are rather pathological and misleading. 

We introduce now an alternative approach to deter- 
mine the parameters for which the system shows col- 
lective firing. Notice that r(t) is time independent in 
regimes (i) (units at rest) and (iii) (desynchronized state) 
but not in (ii). Following previous analysis in the limit 
of infinite units [12[, we write r(t) = pe'* and, assum- 
ing that r(t) is time independent, derive the following 
equations for the global amplitude p and phase >F: 



psin^ = u- um/1 -g{uj)duj = fi{b), (10) 



Kp 2 + pcos^=l \Jb 2 -u 2 g(uj)duj = f 2 (b) 1 (11) 
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with b = \J\ + K 2 p 2 + 2Kpcos*$>. It is possible to ob- 
tain a closed self-consistent relation for b: 



K(f? + fl 



bfr±yjfj + (l-b 2 )f? 



= h±(b). (12) 



For parameters corresponding to regime (ii) Eq. (fT2"|) 
has only one solution corresponding to an unstable fix 
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FIG. 3: Same as Fig[2]for Lorentzian L\(uj) (upper and middle 
panels) and L 4,(01) (upper and lower ones) distributions. Note 
that for the Lorentzian distribution L\(u) the Ott-Antonsen 
method predicts that C, is zero for all values of the diversity A 
while the numerical simulations for finite system size N show 
a rounded-off transition. 



point. When the branches h+(b) and h- (b) meet two new 
solutions arise. This signals a SNIC bifurcation where 
a stable fix point and a saddle are created. Imposing 
h-(b) = h+(b) = b, one can determine the line a (w) 
at which the transition takes place. Noticeably, on the 
SNIC one obtains 'F = 7r/2 independent of the form of 
the distribution g(uj) or the values of the parameters. As 
shown in Fig. 2]the results of this approach coincide with 
those of the Ott-Antonsen method. 

For the Gaussian distribution, the SNIC line that lim- 
its the collective firing regime starts at uJ = 1, a = 0, with 
a negative slope. Therefore for ZU < 1, as a increases one 
finds first a stable steady state for r (regime (i)), then 
crosses the SNIC lower boundary entering in regime (ii) 
and finally crosses the reentrant upper boundary enter- 
ing in regime (iii) where a stable steady state is present 
again. Region (ii) moves upwards and broadens increas- 
ing the coupling K . The overall phase diagram still has a 
Takens-Bogdanov point as in Fig. 1. Decreasing K this 
point moves down and to the left. For low values of K , 
for example K = 1, the upper boundary of the collective 
firing region is limited by the Hopf rather than by the 
reentrant SNIC, nevertheless the diversity-induced tran- 
sition is still present. This reentrance behavior is generic 
for all distributions with well defined moments. 

The Lorentzian distribution behaves differently. The 
SNIC starts at uj = 1,A = with positive slope. Al- 
though apparently this is a small quantitative difference 
and the phase diagrams may be topologically equivalent, 



there are important qualitative consequences and, in par- 
ticular, no diversity-induced transition to collective firing 



14] 



exists 

For distributions which decay fast enough we can ob- 
tain an analytical approximation for the line <r(ui) at 
which the transition to collective firing appears. We as- 
sume that g(u>) ~ for \u\ > b. As b > \Kp — 1|, this 
approximation turns out to be appropriate for large val- 
ues of K. Using vp = 7r/2 one gets p = u; from Eq. 
(|10|). Inserting this in ([TTjl and expanding the integrand 
to second order in oj/b we obtain 



f(2K 2 - 2K\J\ + K 2 ZJ 2 - 1) + 2. (13) 
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From Eq. (I13[) one can obtain da /did. At the origin 
of the SNIC line (u = 1), the derivative takes a nega- 
tive value for any K > and tends to — oo as — K 3 for 
K — > oo. Notice that this is independent of the form of 
g(u) provided the tails decay fast enough. Higher order 
corrections do depend on the specific form g(u>). For the 
Gaussian the next order gives: 
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^JuJ 2 [6uJ 2 - 24K{1 + K 2 uJ 2 )i] + 28(1 + K 2 uJ 2 ) 2 



This approximation, plotted as dashed lines in Fig. 21 
describes well the location of the transition to firing be- 
tween regimes (i) and (ii) for the case of a Gaussian dis- 
tribution. 
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FIG. 4: (Partial) Bifurcation diagram for Gaussian (left) and 
Lorentzian (right) distributions. For the Gaussian, circles cor- 
respond to the Ott-Antonsen method solved numerically as in 
Fig. ([2]), solid lines to the numerical solution of Eas. (|10lTl2)) 
and dashed lines to Eq. (|14p . For the Lorentzian distribution 
the lines correspond to the Ott-Antonsen method. 

In summary, globally coupled active rotators display 
a diversity induced transition to a state of collective 
firing. We have shown that this is a genuine transi- 
tion that is found for any disorder distribution with 



well defined moments. Curiously enough, the transition 
is not present for a Lorentzian distribution, for which 
the first moment integral is only defined as a principal 
value. We have also found that these non-generic results 
given by the Lorentzian distribution also occur in an- 
other well-known excitable system, an ensemble of cou- 
pled FitzHugh-Nagumo 13] units, for which the reen- 
trant diversity-induced transition is present as well for 
distributions such as Gaussian or uniform, but not for 
the Lorentzian one. 

While the Lorentzian distribution seems to be suit- 
able to study generic properties of the original Kuramoto 
model, it is not appropriate for excitable systems and be- 
yond those may not be adequate for other systems. This 
is a clear warning about its use in analytical approxima- 
tions intended to draw conclusions on the generic prop- 
erties of coupled oscillators. 
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For To > 1 there is a small region limited by the homo- 
clinic and the Hopf bifurcations where there is bistability 
between a steady state and one with collective motion. 
However in the collective motion state, r moves slightly 
around a fixed value rather than performing collective 
firing, and the basin of attraction of this state is very 
small, so this is indeed different from the regime we are 
considering. 



